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This paper proposes the application of a novel metaphor-free population 
optimization based on the mathematics of the Runge Kutta method (RUN) 
for parameter extraction of a double-diode model of the unknown solar cell 
and photovoltaic (PV) module parameters. The RUN optimizer is employed 
to determine the seven unknown parameters of the two-diode model. Fitting 
the experimental data is the main objective of the extracted unknown 
parameters to develop a generic PV model. Consequently, the root means 
squared error (RMSE) between the measured and estimated data is 
considered as the primary objective function. The suggested objective 
function achieves the closeness degree between the estimated and 
experimental data. For getting the generic model, applications of the 
proposed RUN are carried out on two different commercial PV cells. To 
assess the proposed algorithm, a comprehensive comparison study is 
employed and compared with several well-matured optimization algorithms 
reported in the literature. Numerical simulations prove the high precision 
and fast response of the proposed RUN algorithm for solving multiple PV 
models. Added to that, the RUN can be considered as a good alternative 
optimization method for solving power systems optimization problems. 
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1. INTRODUCTION 


The increased demand for electrical energy while decreasing CO2 emission associated with 
electricity generation requires green renewable energy sources (GRERs). As a source of green energy, 
photovoltaic (PV) can generate electricity from sunlight using semiconductor materials. Globally, there is a 
massive increase in solar photovoltaic capacities in comparison with those of wind energy as per the REN21 
latest published report in 2021 [1]. The Egyptian government pays a concentration towards the exploitation 
of available GRERs and PV especially. One of the most significant constructed PV projects in the world was 
at Benban, Aswan, Egypt, which has an installed capacity of 1.8 GW [2]. 

Due to this high potential spread, accurate modeling of PVs has a high impact on their performance 
under varied environmental conditions and shading. To express the non-linearity relationship between 
current-voltage (I-V) and power-voltage (P-V) of the PV cell, equivalent circuit models will be required. As 
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per the literature, the single and double diode equivalent circuit models are used in parameter extraction of 
PV cells, which can be reviewed in [3]. However, its simplicity and reduced number of estimated parameters 
of the single diode equivalent circuit model, its precision diminishes at low irradiation levels with 
temperature variation [4], [5]. The estimated parameters in the single diode model are the photogenerated 
current, the reverse saturation currents, the series and shunt resistances, and the ideality factor. The double 
diode equivalent circuit model is composed of seven parameters, five of them the same as five parameters of 
the single diode in addition to two parameters, which are the reverse saturation current and the ideality factor 
of the additional diode. 

As per the literature review, the parameter extraction techniques of the PV equivalent circuit model 
can be classified into three categories, namely analytical, numerical (deterministic and stochastic), and hybrid 
techniques [6], [7]. The analytical methods utilize the power-voltage and current-voltage data curves with the 
aids of selected data points from the manufacturer datasheets (the open-circuit voltage V_oc, the maximum 
output power P_m, the short circuit current I_sc, the maximum output current I_m, the maximum output 
voltage V_m) to constitute the mathematical parameter estimation problem [8]-[13]. These methods require 
some mathematical approximations (simplification) to reduce the number of extracted parameters. Although 
the used approximation provides ease of implementation and less computational effort, however, it has a 
significant impact on the solution accuracy. In study [14], it was approved that these approaches are less 
accurate than numerical approaches. To overcome modeling imprecision that may arise if there are inaccurate 
selected data from datasheets, the measured data of I-V are the same as used in numerical (iterative) 
techniques. The parameter extraction techniques can be categorized into two categories: deterministic and 
stochastic (heuristic and meta-heuristic) techniques. The deterministic techniques such as Newton-Raphson 
[15], linear identification [16], or the Levenberg-Marquardt (LM) algorithm [17] can be used for parameter 
extraction of PV equivalent circuit parameters. The main drawback of these conventional methods is local 
optimum trapping due to its sensitivity to the initial solution. As an alternative to the deterministic 
techniques, the naturally inspired (meta-heuristic) are extensively used in the last decade. In literatures, there 
are a lot of these methods such as, biogeography-based heterogeneous cuckoo search (BHCS) [18], pattern 
search (PS) [19], firefly algorithm (FA) [20], ant lion optimizer (ALO) [21], Jaya algorithm [22], salp swarm 
algorithm (SSA) [23], elephant herd optimizer [24], enhanced sine cosine algorithm (ISCA) [25], hybridized 
interior search algorithm (HISA) [26], an artificial bee colony-differential evolution (ABC-DE) [27], 
improved adaptive Nelder-Mead simplex (NMS) hybridized with the artificial bee colony (ABC) 
metaheuristic, algorithm of hybrid adaptive and Nelder-Mead simplex (EHA-NMS) [28], mutative-scale 
parallel chaos optimization algorithm (MPCOA) [29], classified perturbation mutation based particle swarm 
optimization (CPMPSO) [30], heterogeneous comprehensive learning particle swarm optimizer (HCLPSO) 
[31], and improved shuffled complex evolution (ISCE) [32]. A forensic based optimization algorithm was 
developed in [33] for finding the optimal parameters of solar cell modules. Another optimizer called 
turbulent flow of water optimizer in [34] was developed to optimize the parameters of three solar cell models. 
An assessment study based on the elephant herd optimization, which developed with different versions [35], 
is compared with closed loop particle swarm optimizer in [36]. The An interval branch and bound global 
optimization algorithm (IBBGO) that is referred to, interval branch/bound global optimizer, was integrated to 
find the optimal parameters of three PV models [37]. 

All these meta-heuristic optimization algorithms aim at minimizing the objective/cost function while 
balancing between exploration and exploitation phases. At the same time, meta-heuristic optimization 
algorithms cannot provide robust search capability towards optimal solutions [38]. To overcome these 
challenges, a novel metaphor-free population optimization based on the mathematics of the Runge Kutta 
method (RUN) is presented [39]. RUN optimizer balances between exploration and exploitation phases in 
dynamic behavior. Moreover, the RUN optimization algorithm has a competitive convergence speed as well 
as an enhanced solution quality to avoid local optimal solutions. 

It is cleared that many meta-heuristic optimization algorithms are used for parameter extraction of 
the equivalent circuit of the PV model. Nevertheless, most of these algorithms have several control 
parameters that require tuning to achieve better performance for an optimization problem. The RUN 
algorithm has fewer control parameters due to its specific property which depends on the main logic of the 
Runge Kutta technique. This advantage motivates the authors of the presented article to use it in parameter 
estimation of the double diode model of PV cell/module. The main contribution of this research can be 
summarized as: i) develop RUN optimization algorithm for PV parameter extraction based on double diode 
model, ii) assessment of the extracted parameter using RUN with the recent optimization algorithms, iii) 
recommending the best method and best equivalent circuit model to be used for PV cell/module, and iv) the 
proposed RUN has the best convergence rates compared with the competitive methods. 

This paper will be organized as follows: After this section, problem formulation will be introduced 
in section 2. Section 3 will develop the RUN optimization algorithm for use in parameter extraction of a PV 
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cell/module. Application of the proposed algorithm to different commercial cell/module and test results will 
be emphasized in section 4. Section 5 will conclude the results drawn from this research. 


2. PROBLEM FORMULATION 
2.1. Electrical models 

Mathematical modeling of PV cell/module equivalent circuits will be illustrated. The equivalent 
circuit models emulate the non-linear (I-V) and (P-V) relationships of the PV cell. This model is considered 
as the simplest equivalent circuit model. It has five parameters which are the photogenerated current Ipn, 
series and shunt resistances R,, R,,, the ideality factor n and the reverse saturation current Ip, as shown in 
Figure 1(a). The load current I of solar cells can be represented by (1): 


Taj I | (- +IRs 1)] (V + IR;) 
= _ ~~ ex — — — —— 
OTRA Ron (1) 

In this model, the recombination of generated charge carriers was neglected. Moreover, its precision 
diminishes at low irradiation levels with temperature variation [4]—[6]. The recombination effect was 
modeled in the double diode model that is presented in Figure 1(b), the load current J of solar cells can be 
represented by (2): 


taj I | (- +IR; 1)| I | (- + IR, 1)| (V + IR;) 
= — Ips_1 |ex — — Ięgs-2 | ex — -> 
ph — trRs-1 |EXP nV; RS-2 p Ro, (2) 


nV; 


ph 


(a) 


Figure 1. Diode equivalent circuit models, (a) single and (b) double 


2.2. Mathematical representation of parameter extraction problem 

The parameter extraction problem is considered an optimization problem. The main objective of the 
estimation problem is to minimize the variance between the experimental data and simulated results so that 
the optimal values of these unknown model parameters can be extracted. The objective function is defined as 
the overall root mean square error (RMSE). For N-measurements, the objective function is formulated as [38] 
as in (3): 


1 N 
RMSE (x) = =Y Velox) — I]? 
=] (3) 


where, x = Hoi Irs-1: I ps—2, Rss Rsn: M1) nz} and the f (Vp, Ig, X) is used for current calculation from (2). 


3. RUNGE KUTTA OPTIMIZER (RUN) 

The Runge Kutta-based optimization (RUN) is a novel optimization algorithm proposed in [39]. The 
RUN is developed based on the mathematics of the Runge Kutta method. It is inspired by the logic of slope 
variations computed by the Runge Kutta method as a promising and logical searching mechanism for global 
optimization. In RUN, three distinctive phases are simulated. These phases are active exploration and 
exploitation phases for exploring the promising regions in the feature space and constructive toward the 
global best solution. Moreover, in the last phase, an enhanced solution quality mechanism is employed to 
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avoid the local optimal solutions and increase convergence speed. Figure 2 demonstrates the phases of the 
proposed RUN algorithm. The phases of the suggested RUN optimizer are discussed in detail in the 
following sub-sections. 


Algorithm 1. The pseudo-code of RUN 
Stage 1. Initialization 
Initialize a, b 
Generate the RUN population X, (n = 1, 2, ..., N) 
Calculate the objective function of each member of the population 
Determine the solutions xy, Xp, and Xpest 
Stage 2. RUN operators 
for i= 1: Maxi 
forn=1:N 
for/=1:D 
Calculate position x, 41, using Eq. 4 
end for 
Enhance the solution quality 
if rand < 0.5 
Calculate position Xnew2 using Eq. 7 
if f (xn) < f Œnew2) 
if rand< w 
Calculate position Xnew3 using Eq. 9 
end 
end 
end 
Update positions xw and xp 
end for 
Update position Xpest 
i=i+1 
end 
Stage 3. return Xpest 


Figure 2. Pseudo-code of the proposed algorithm (RUN) 


3.1. Updating solutions stage 
In the updating solutions stage, the RUN uses a search mechanism (SM) based on the runge kutta 
method to update the position of the current solution at each iteration, which is defined as (4): 


if rand < 0.5 
(Exploration phase) 
Xn+1 = (xe + r. SF. g.Xc) + SF.SM + u. randn.(Xm — xc) 
else 
(Exploitation phase) 
Insi = Om TSF. g-m) + SF.SM + p.randn. (Xi = Xr) 
end (4) 


where, r is an integer number, which is 1 or-1. g is a random number in the range [0, 2]. SF is an adaptive 
factor. Where, u is a random number. The formula of SM is defined in [34]. The formula of SF is as (5): 


SF = 2.(0.5 — rand)xf (5) 


where, f = a X ex p (-» x rand x ()). Maxi stands for the largest number of iterations. The formula 


of x: and X» are as (6): 


Xe =QOQX%y,+A-@Q) XxX 
Xm = P X Xbest t (1 a ~) X Xipest (6) 


where, @ is a random number in the range of (0,1). xpesr is the best-so-far solution. xipesr is the best position 
obtained at each iteration. 


3.2. Enhanced solution quality stage 

In the RUN algorithm, enhanced solution quality (ESQ) is employed to increase the quality of 
solutions and to avoid local optima in each iteration. The following scheme is executed to create the solution 
(Xnew2) by using the ESQ (7), (8): 
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if rand < 0.5 

ifw<l1 

Xnew2 = Xnew1 + r.WA Snowe 5 Xiva) + randn| 

else 

Xnew2 7 (Xnew1 a Xavg) + r.W.|(u.Xnew1 a Xg) + randn| 

end (7) 
= err eee = Mra Ara rs 

w = rand (0, 2).exp ( c ()) avg 5 3 

Xnew1 = Ê X Xavg + (1 — P) X Xpest (8) 


where, p is a random number in the range of [0, 1]. c is a random number, which is equal to 5X rand in this 
study. r is an integer number, which is 1, 0, or-1. Xpest is the best solution explored so far. The solution 
calculated in this part (Xnew2) may not have better fitness than that of the current solution (i.e., f (Xnew2) > 
f(xn)). To have another chance for creating a good solution, another new solution (Xnew3) is generated, 
which is defined as in (9), 


if rand < w 
Xnew3 = (Xnew2 ~ rand.Xnew2) + SF.(rand.xpx + (v.xp a Xnew2)) 
end (9) 


where v is a random number with a value of 2 x rand. 


4. SIMULATION AND RESULT 

To assess the proposed algorithm, a comprehensive study is employed compared with several 
previous optimization algorithms. Numerical simulations on different commercial PV cells/modules will be 
illustrated below. To evaluate the performance of the proposed algorithm, A set of standard data for a 
commercial silicon solar cell (made by R.T.C. company from France) with a diameter of 57 mm, at a 
temperature of 33 °C, and 1 sun (1,000 W/m?) [40]. The PhotowattPWP201 PV commercial module with 36 
polycrystalline silicon cells connected in series, operating under an irradiance of 1,000 W/m? and temperature 
of 45 °C [41]. Table 1 summarizes the datasheets of the selected commercial PV cell/module to be tested in 
this work. The boundary constraints of the extracted parameters are given in Table 2. The proposed RUN- 
based model and the other comparative algorithms are executed by the authors via the matrix laboratory 
(MATLAB) 2017a platform using an Intel ® core TM i5-7200U CPU, 2.50 GHz, 8 GB RAM Laptop. 
Accuracy examination of the proposed optimization RUN algorithm for parameters identification will be 
additionally accomplished by the current calculation based on the values estimated for the two models 
considered for comparison with that taken from the experimental measurements. The error concerning each 
of the measured values was evaluated by relative error (RE) and individual absolute error (IAE), calculated 
as given in (10) and (11), respectively. 


RE = (Imeasured a lictimated)/ Imeasured (10) 


IAE = lIMeasurea > Testimatea| (11) 


Table 1. Datasheets of selected cell/module for study 


Company Cell type Voce LV] Isc [A] Vn [V] In [A] Pm [W] Reference temp. [C] _ # of cells / module 
France solar NA 0.5728 0.76 0.45 0.691 11.315 33 1 
Photowatt Polycrystalline 16.778 1.03 12.60 0.898 0.311 45 36 


Table 2. Boundary constraints of the extracted parameters 


Model Parameters 
France solar Ipn [A] Irs1,2[4A] Rs [Q] Rsn [Q] Ni» 
Min 0 0 0 0 1 
Max 1 1 0.5 100 2 
Photowatt 201 Min 0 0 0 0 1 
Max 2 10 0.5 3000 2 
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4.1. RTC. france solar cell 

For RTC solar cell model, the optimal parameters attained using different optimization techniques 
are presented in Table 3. The solar cell parameters estimated using the RUN algorithm have the lowest 
RMSE value. Figure 3 demonstrates that the estimated I-V and P-V characteristics using RUN are identical 
compared to the measurement as well as the datasheet. Table 4 spots the light on the coincidence between 
both the measured and estimated values for voltage and current. The results of MPP are highlighted in green 
in Table 4 and illustrated in Figure 3. As can be concluded from Figure 4, RUN has the lowest RMSE value 
(0.0009829) in comparison with sunflower optimization (SFO), which has the highest value (0.002). On the 
other hand, RUN has the lowest execution time of (9.95 s) in comparison with SFO, which has the highest 
execution time of (44.402 s). 


Table 3. Optimum parameter settings using different optimization techniques for RTC France solar cell 


RCGA CSA SSA PSO SFO GW-CS RUN 
Ipn (A) 0.7606 0.7608 0.7604 0.76077 0.76275 0.76002 0.76077 
Tps—1 (HA) 0.3743 0.6489 0.3165 0.6680 0.2452 0.2272 0.22266 
Tps—2 (HA) 0.0731 0.2199 0.1728 0.2389 0.8324 0.2553 0.67737 
R; (Q) 0.0357 0.0367 0.0369 0.0366 0.0339 0.0359 0.03674 
Rsa (Q) 60.823 54.946 59.7711 55.691 51.165 67.0574 55.4973 
ny 1.4968 1.9454 1.69629 2.000 1.4768 1.4632 1.45006 
Nz 1.9607 1.4493 1.43602 1.4559 1.7895 1.6630 1.96604 
RMSE 0.0010 0.000985 0.00103 0.00098 0.002 0.00110 0.0009829 
Duration [s] 13.655 11.524 10.1512 10.567 44.402 14.131 9.95 
0.8 T T T T T T T T T T 04 
0.75 = 
0.7 7 
0,65 Pe | I ES J 
0.6 = a 
0.55 I mee. 
0.5 l 4 
L045 l TS 
5 0.4 I 40.2 g 
0.35 I J Å 
0.3 l | 
0.25 l = 
0.2 I 40.1 
0.15 I a 
se = Estimated Current © Measured Current O Measured Power l | 
0.05 v a] 
0 | i i L 1 | | i 1 y MpP | fi o 
0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 
Voltage (V) 


Figure 3. The coincidence between measured and estimated results in V-I and V-P curves of RTC-solar cell 


Table 4. The relative error for each measurement using a RUN-based two diode model 


Measured Measured Estimated IAE RE Measured Measured Estimate IAE RE 

Voltage current current Voltage current d current 

-0.2057 0.764 0.7640 0.0000 0.0000 0.4373 0.728 0.7399 0.0014 0.0019 
-0.1291 0.762 0.7626 0.0006 0.0008 0.4590 0.7065 0.7271 0.0009 0.0012 
-0.0588 0.7605 0.7613 0.0008 0.0011 0.4137 0.6755 0.7066 0.0001 0.0001 
0.0057 0.7605 0.7602 0.0003 0.0004 0.4784 0.632 0.6748 0.0007 0.0010 
0.0646 0.7600 0.7591 0.0009 0.0012 0.4960 0.573 0.6302 0.0018 0.0029 
0.1185 0.7590 0.7581 0.0009 0.0012 0.5119 0.499 0.5711 0.0019 0.0033 
0.1678 0.7570 0.7572 0.0002 0.0002 0.5265 0.413 0.4982 0.0008 0.0017 
0.2132 0.7570 0.7562 0.0008 0.0010 0.5398 0.3165 0.4118 0.0012 0.0030 
0.2545 0.7555 0.7552 0.0003 0.0004 0.5521 0.212 0.3150 0.0015 0.0046 
0.2924 0.7540 0.7537 0.0003 0.0004 0.5633 0.1035 0.2095 0.0025 0.0119 
0.3269 0.7505 0.7514 0.0009 0.0012 0.5736 -0.01 0.0997 0.0038 0.0368 
0.3585 0.7465 0.7473 0.0008 0.0010 0.5833 -0.123 -0.0126 0.0026 -0.2618 
0.3873 0.7385 0.7640 0.0000 0.0000 0.5900 -0.21 -0.1280 0.0050 -0.0408 
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Figure 4. The objective function of the competitive estimation algorithms 


4.2. Photo watt PWP201 

In the Photowatt PWP201 module, from the obtained results in Tables 5 and 6, the RUN has the 
lowest RMSE value (0.003139) among the competitive optimization algorithms. The maximum power point 
(MPP) is highlighted in green in Table 6 and illustrated clearly in Figure 5. The RMSE for all optimization 
techniques are plotted in Figure 6. Generally speaking, RUN has the lowest RMSE value (0.003139) 
compared to all the applied optimization techniques. The crow search algorithm (CSA) comes in the first 
class in terms of execution time with (10.2910 s) in comparison with the SFO with (46.7350 s), which has the 
highest value of elapsed time. 


Iterations 


998 


Iterations 


998.1 


1000 


Table 5. Optimum parameter settings using competitive optimization algorithms for Photowatt-PWP201 


RCGA CSA SSA PSO SFO GW-CS RUN 

Ipn (A) 1.0252 1.0258 1.02867 1.0258 1.01883 1.02574 1.0270 

I ps1 (HA) 5.1463 9.8979 7.59250 3.8226 0.2769 0.2527 1.4074 

Irs-2(HA) 2.4599 3.8944 7.53697 5.0705 0.9499 6.01257 5.5321 

R; (Q) 0.0321 0.0322 0.02864 0.03175 0.03724 0.03191 0.0321 
Rn(Q) 1000 1001 629.437 1000 351.092 248.962 100.3905 

ny 1.3982 1.877 1.48011 1.7271 1.32756 1.89869 1.9920 

n, 1.767 1.3723 1.61563 1.4007 1.23380 1.41230 1.4034 
RMSE 0.0035 0.0035 0.00491 0.00349 0.00734 0.00347 0.003139 

Duration[s] 12.0730 10.2910 19.5530 10.4660 46.7350 11.4380 11.25 


Note: RCGA= real coded genatic algorithm, SSA = introduced in the introduction, PSO is particle swarm optimization, 
GW-CS = gray wolf cuckoo search algorithm 


Table 6. The relative error for each measurement using a RUN-based two diode model 


Measured Measured 


Estimated 


Measured 


Measured 


Estimated 


IAE RE IAE RE 

voltage current current Voltage current current 

1.9426 1.0345 1.0272 0.0073 0.0070 12.6490 0.9120 0.9227 0.0028 0.0030 
0.1248 1.0315 1.0266 0.0049 0.0047 13.1231 0.8725 0.9117 0.0003 0.0003 
1.8093 1.0300 1.0261 0.0039 0.0038 14.2221 0.7265 0.8715 0.0010 0.0012 
3.3511 1.0260 1.0256 0.0004 0.0004 14.6995 0.6345 0.7264 0.0001 0.0002 
4.7622 1.0220 1.0249 0.0029 0.0029 15.1346 0.5345 0.6352 0.0007 0.0010 
6.0538 1.0180 1.0239 0.0059 0.0058 15.5311 0.4275 0.5349 0.0004 0.0007 
7.2364 1.0155 1.0222 0.0067 0.0066 15.8929 0.3185 0.4285 0.0010 0.0024 
8.3189 1.0140 1.0189 0.0049 0.0049 16.2229 0.2085 0.3188 0.0003 0.0010 
9.3097 1.0100 1.0131 0.0031 0.0031 16.5241 0.1010 0.2083 0.0002 0.0009 
10.2163 1.0035 1.0029 0.0006 0.0006 16.7987 -0.0080 0.0989 0.0021 0.0204 
11.0449 0.9880 0.9862 0.0018 0.0018 17.0499 -0.1110 -0.0076 0.0004 -0.0472 
11.8018 0.9630 0.9603 0.0027 0.0028 17.2793 -0.2090 -0.1106 0.0004 -0.0034 
12.4929 0.9255 1.0272 0.0073 0.0070 17.4885 -0.3030 -0.2091 0.0001 -0.0005 
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Figure 5. The coincidence between measured and estimated results in V-I and V-P curves of PWP201 
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0 100 200 300 400 500 600 700 800 900 1000 
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Figure 6. The convergence rates of objective function decline using competitive algorithms 


5. CONCLUSION 

This manuscript presents a new physical implementation of modern inspired optimization 
algorithms called RUN for PV cells/modules parameter estimation. To validate and verify the main findings, 
the proposed algorithm was implemented and tested on two commercial cells/modules. Moreover, a 
comprehension comparison of the latest meta-heuristic optimization algorithms was illustrated. The attained 
theoretical and experimental results are coincident, which proves the superiority of RUN in the field of 
parameter estimation for PV cells/modules. The results signify the effectiveness and reliability of the 
proposed RUN in estimating the accurate double diode model of two practical PV cells/modules. The RUN 
realizes steady convergence rates than other competitive algorithms. Finally, the main findings of this article 
will pave the way for the authors as well as the researchers to evaluate the impact of the estimated parameters 
in emulating PV equivalent circuits under various real operating conditions. 
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